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Abstract. The distribution and occurrence of a species in its habitat are inevitably linked to its ecology. To successfully 
monitor and protect species, it is important to investigate which species-specific factors influence its interactions with the 
environment. In this study, we focus on patterns in habitat use of the Sand Lizard (Lacerta agilis). Differences in seasonal 
as well as sex and size class-dependent habitat use have been reported from the edges of this species’ range. To verify such 
trends in the core area of its distribution, we analyzed the habitat factors weather, microclimate, microhabitat structures, 
and time dependence, which may have an impact on the use of space of the Sand Lizard. Using generalized linear models, 
hypervolumes, density estimations, and Chi-squared tests, we found that the movement patterns of individuals can neither 
be described by time differences, climatic conditions, or habitat composition, nor do they show habitat- or weather-related 
differences of movement between sexes or size. Here we demonstrate that in the case of a population from the core of this 
species’ distribution area in the Dellbriicker Heide (Germany), habitat use is solely influenced to a low degree by differ- 
ences related to the ontogeny of Sand Lizards and does not depend on any of the other evaluated factors. These results are 
in enormous contrast to findings in populations from the periphery of their distribution, i.e., the United Kingdom, Latvia, 
Romania, Bulgaria, and the Pyrenees. This implies that seasonal habitat shifts are more extreme at the edges of the range 


of L. agilis, serving to compensate deteriorating habitat conditions in the periphery. 
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Introduction 


Geographic distributions of species are determined by the 
ecological niches they occupy, which in turn are based on 
abiotic factors, biotic interactions, geographic accessibil- 
ity, and characteristics of the individuals themselves (So- 
BERON & PETERSON 2005). The fundamental niche of a 
species comprises all parameters and resources, which are 
necessary for an indefinite maintenance of viable and re- 
producing populations and so determine their potential 
distribution. The realized niche of a species is shaped by 
the intersection of biological interactions (e.g., predation 
pressure, competition, resource availability), environmen- 
tal factors such as geographic accessibility, and the physi- 
ological potential of the species itself (HUTCHINSON 1957, 
GRINNELL 1917, ELTON 2001). 

KUHNELT suggested in 1943 that habitat choice parame- 
ters of species with a particularly wide distributional range 
can vary enormously depending on the geographic loca- 
tion. The diversity of habitat types that can be colonized 
by the species should be larger in the core area, where con- 
ditions are optimal, and the species here becomes eury- 


biontic, which suggests that individuals can tolerate a 
broad range of a given habitat factor (SCHAEFER 2012). At 
the edges (periphery) of large distribution ranges, condi- 
tions typically become more compromised due to a re- 
duced availability of suitable microhabitats and hence force 
the lizards to be more selective. This is the reason why the 
species regionally becomes stenobiontic, suggesting that 
it can tolerate only a narrow range of a given habitat fac- 
tor (KUHNELT 1943, BOHME 1978, SCHAEFER 2012, BOHME 
& RODDER 2014). Next to overall habitat choice across a 
species’ range, seasonal habitat shifts within single popula- 
tions can be linked to the same principle. 

The Sand Lizard, Lacerta agilis LINNAEUS, 1758, has one 
of the most extensive distribution ranges of all the world’s 
reptile species, but has, in some areas, suffered heavy de- 
clines especially in the northwestern portions of its distri- 
bution (EDGAR & BIRD 2006). Main threats to this species 
are habitat fragmentation and destruction of microhabitat 
structures due to human activities (EDGAR & BIRD 2006, 
AGHASYAN et al. 2021). Additionally, some populations at 
the edges of its distribution are quite fragmented and iso- 
lated, like in the north and northwest of the United King- 
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dom where populations have declined due to habitat loss 
(HOUSE & SPELLERBERG 1983, WOODFINE et al. 2017). Some 
isolated populations occur in the mountains of Greece, 
Bulgaria, and in the Pyrenees at altitudes of up to at least 
2500 m a.s.l. (Bulgaria). They can be found in these re- 
gions only in highly specific habitats with suitable (micro-) 
climatic conditions (BOHME 1978, EDGAR & BIRD 2006, 
AGHASYAN et al. 2021). 

At its northern range limits, the Sand Lizard is re- 
stricted to habitats with warm and dry climates suitable 
for thermoregulation, wherein this species tends to use 
colder and more humid conditions than it does on aver- 
age at its southern range limits. Therefore, we expected 
these lizards to choose their habitats accordingly, based on 
available habitat types. In Latvia, L. agilis prefers vegeta- 
tion composed of sparse low pines (Pinus sylvestris) with a 
high degree of herbal ground cover dominated by grasses 
and heath. The lizards generally prefer dry areas and avoid 
herbal growth with tall broad leaves (CErRANS 2007). In 
Romania, researchers found age differences in habitat se- 
lection, with adult and yearling lizards tending to use the 
compact vegetation at the bases of bushes whereas hatch- 
lings were found using lower vegetation and open patches 
(Nemes et al. 2006). Sexual and ontogenetic differences in 
microhabitat use and phenology have previously been re- 
ported and are related to reproductive cycles (SAINT-GI- 
RONS 1976, NULAND & STRIJBOSCH 1981). In the Pyrenees, 
it was found that Sand Lizards exhibited differences in 
habitat use and activity, depending on sex and time of the 
season, with meadows being most commonly used, except 
during breeding season, when these lizards preferred bush- 
es over meadows. The authors also found that juveniles 
rarely used bushes and preferred meadows that had been 
cropped short (Amar et al. 2003). In Bulgaria, researchers 
suggested that the height and density of grass seems to be 
a key microhabitat feature that spaces the sand lizard’s dis- 
tribution in the subspecies L. a. bosnica (PRIETO-RAMIREZ 
et al. 2018). Furthermore, differences in microhabitat selec- 
tion were detected between males and females of the sub- 
species L. a. chersonensis (GROZDANOV et al. 2014). In the 
United Kingdom, which is part of the northwestern edge 
of the Sand Lizards’ range, this species inhabits exclusively 
sandy and dry environments like coastal dunes or heath- 
lands and prefers more open habitats (HOUSE & SPELLER- 
BERG 1983, WOODFINE et al. 2017). JACKSON (1978) already 
hypothesized that climate may have an influence on Sand 
Lizard distribution in Britain. It was suggested that a low 
level of sunshine in May was one of the reasons for the 
decline of the Sand Lizard population on the Merseyside 
coast in northwestern England in the 1960’s (JACKSON 
1978). These findings from the edges of the Sand Lizard’s 
range show up extreme variations in habitat choice, mak- 
ing peripheral populations seem to compensate for overall 
poorer conditions by responding more extremely to differ- 
ent parameters of their habitat. Studies performed in these 
peripheral situations already provide important informa- 
tion on habitat selection and habitat use by this species, 
including that they can vary drastically between individual 


locations. However, to our knowledge, no study examin- 
ing factors influencing spatio-temporal habitat selection 
in the core distribution area of the Sand Lizard has as yet 
been conducted. Thus, based on the available data from 
the periphery we hypothesize: As imposed by the principle 
of KUHNELT (1943), we expect spatio-temporal habitat use 
in the core distribution area of the species to be less pro- 
nounced compared to the range edges. To test this hypoth- 
esis, we monitored conditions in the distributional core 
area of L. agilis, focusing on the nature reserve Dellbriicker 
Heide, North Rhine-Westphalia, Germany. We here test 
if distances to habitat structures and the microhabitat use 
of Sand Lizards are subject to seasonal or daily differences 
due to changes in weather conditions. Furthermore, we as- 
sess if these parameters as well as microclimatic conditions 
differ depending on the sex and/or size of individuals. Our 
results from the core area are compared to data from the 
Sand Lizard’s range edges. 


Materials and methods 
Study area 


The nature reserve Dellbriicker Heide (geographic location 
in WGS 84: 50.98187° N; 7.05805° E, altitude: 37-59 m a.s.l) 
with an area of only 40 ha, is a very small part of the Bergi- 
sche Heideterrasse (KULADIG 2016). The Bergische Heide- 
terrasse has a total length of 80 km, ranging from Duisburg 
to Siegburg, and was formed from the mountain fuselage 
of the Rhenish Slate Mountains by the erosion of older rock 
layers and a deposition of loose more recent layers. De- 
spite its small extent, the Dellbriicker Heide has a remark- 
ably high diversity in flora and fauna composed of sandy 
grass- and heathland areas framed by small forests, which 
provide a perfect mosaic of transition structures. These at- 
tributes are known to be very important for the viability of 
many species (BRACHET et al. 1999) and constitute a highly 
suitable habitat for L. agilis (BLANKE 2021). 


Study design 


To assess the spatial positioning of Sand Lizards in their 
habitat in a standardized and constant manner, visual en- 
counter surveys were conducted using predefined tran- 
sects and random points, covering all available habitat 
types (Fig. 1). Detectability of lizards was not only depend- 
ing on habitat composition, but also on weather conditions 
during our fieldwork, for which reason data collection was 
extremely difficult on days with heavy rain or extreme 
heat. Days with overall poor sampling conditions therefore 
were subsequently avoided by checking the weather fore- 
cast every day prior to a survey (for baseline conditions for 
surveys see Supplementary Fig. S1). Nonetheless we test- 
ed whether sampling days represented the overall weath- 
er conditions in an appropriate manner by using the func- 
tion sm.density.compare from the R package sm (BOWMAN 
& AZZALINI 2018). We created density estimation plots 
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showing the density of all days for weather variables, maxi- 
mum temperature (°C), relative humidity (%), mean wind 
speed (m/s), sunshine, and precipitation duration (hours) 
at ten-minute intervals. Additionally, permutation tests of 
equality were computed to assess their statistical signifi- 
cance (WILCOX 2010), and p-values were analyzed to check 
if sampling days represented random samples. 

Field data were recorded at least 40 times a year (52 
times in 148 days with 66 h of data acquisition in total from 
2 May to 26 September 2018, 40 times in 133 days with 60 h 
of data acquisition in total from 24 April to 3 September 
2019, and 42 times in 97 days with 118 h of data acquisition 
in total from 7 June to 11 September 2020). Due to natural 
changes of the environment, such as vegetation becoming 
impenetrable, landslides, or other barriers, some small ad- 
justments to the borders of transects and numbers of ran- 
dom points were necessary. 

All transects and random points were visited equally 
alternating at different time intervals, at least once in the 


morning from approximately 9:00 to 12:00 h, at midday 
from approximately 12:00 to 15:00 h, and in the afternoon 
from approximately 15:00 h to 18:00 h, to cover possible 
daily variations in spatial habitat use. Assessments were 
conducted systematically by following a predefined and 
memorized path based on prominent landmarks like indi- 
vidual shrubs and trees to avoid pseudo-replications. Ran- 
dom points were searched within a radius of about 20 m 
depending on accessibility. The following data were col- 
lected for every Sand Lizard observation: 

First, the GPS coordinates of the sampling spot was tak- 
en in decimal degrees (WGS 84) with a maximum of 5 m 
uncertainty in radius. Furthermore, the date and exact 
sampling time were noted. Microclimate data included air 
temperature (in °C to the nearest 0.1°C + 3.0%) and relative 
humidity (in % rH to the nearest 0.1% + 5.0%) (multimeter 
ELV VA19), approximately 70 cm above the sampling spot 
(for summary statistics of the measurements see also Sup- 
plementary Table S2), cloud cover in three categories, o = 


Figure 1. Aerial imagery of the study area Dellbriicker Heide. Red lines indicate transect locations. Red abbreviations identify the 
different transects after compass directions. White dots indicate random point locations with abbreviations as individual identifiers 
of each point. 
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no clouds, 1 = some clouds, 2 = partially cloudy, 3 = closed 
cloud cover, substrate temperature at the sampling spot 
(measured with a laser infrared thermometer ennoLogic 
eT1050D Dual Laser), and brightness of light at the sam- 
pling spot (categorized as sunny [spot fully exposed to the 
sun], half-shaded [shaded with insolated patches or sunny 
with shaded patches], or shaded). Microhabitat structure 
data were assessed within a 5-m radius, including the veg- 
etation and predominant substrate type (sand, soil, gravel, 
or a combination of these). Sex and size classes (male, fe- 
male, subadult) were visually assigned, by rough size esti- 
mations and evaluations of the colour pattern. The category 
‘subadult’ comprises immature juveniles and yearlings that 
still had smaller sizes than adults and had not yet developed 
their typical distinct dorsal pattern and slim head shape of 
females or the green flanks and broader head of males, ei- 
ther of which are visible already from a distance (BISCHOFF 
1984, OLSSON 1988, EDGAR & BIRD 2006). 


Weather data 


Data on the weather conditions within the sampling peri- 
ods were obtained via the German Meteorological Service 
(DWD) from the meteorological station with the ID 2667 
situated at 50°51’ N, 7°09’ E. These included the current air 
temperature (in °C) at 5 cm and 2 m above the ground, 
relative humidity (in %), amount (in mm) and duration of 
precipitation (in hours), solar irradiance (in joule), hours 
of sunshine, maximum, minimum, and mean wind speed 
(in m/s) given for every ten minutes. Here, only the data 
consistent with encounter times of data assessments from 
2 May to 26 September 2018, 24 April to 3 September 2019, 
and 7 June to 11 September 2020, were used for statistical 
analyses. Summary statistics of these weather conditions 
can be found in the Supplementary Table S3, characteriz- 
ing the baseline conditions applied to this study. 


Distance data 


To assess habitat composition, seven classes of land cover 
were defined based on local conditions and our experience 
from the field from the years 2018 and 2019: Water, sand, 
grass (vegetation up to 50 cm in height), scrub (vegetation 
from 50 cm to 5 m in height), trees (vegetation > 5 m in 
height), heath and blackberry. The latter two are treated as 
separate categories due to their being structurally different 
from other scrub and because they must be regarded as im- 
portant structures in Sand Lizard habitats (EDGAR & BIRD 
2006, CLEMENT et al. 2022). 

To create a categorical habitat map, the methods estab- 
lished by CLEMENT et al. (2022) were applied in a modified 
manner. A supervised land cover classification was cre- 
ated based on the latest UTM tiles (2019), available from 
the geoportal of North Rhine-Westphalia (RGB images, 
10 cm pixel resolution). A habitat map was created using 
MaxEnt version 3.4.1. software (PHILLIPS et al. 2017; for 


details, see below), as well as the geographic information 
system QGIS version 3.14.16 (https://qgis.org/de/site/; last 
accessed 15 November 2021), and the statistics program R 
version 4.0.5 (R Core Team 2021). The three spectral bands 
(red, green, blue) of the UTM tile were extracted and used 
as environmental variables. As training records, 100,000 
randomly generated training points were created of which 
as many points as possible were manually assigned to the 
seven habitat classes (water, sand, grass, heath, blackberry, 
scrub, trees). In total, 100 models were computed using a 
bootstrap approach with an 80:20 split for model training 
and testing with standard settings. The resulting probabil- 
ity maps represented the likelihood of each grid cell to rep- 
resent one of the habitat types. Performance of the mod- 
els was evaluated using the Area Under the receiver op- 
erating characteristic Curve (AUC, ,), which evaluates the 
prediction accuracy of the model (Xu et al. 2019). AUC pu 
values above 0.7 indicate a useful discrimination ability of 
the model (PHILLIPs et al. 2017, WALDEN-SCHREINER et al. 
2017, Xu et al. 2019). In a next step, using R (R Core Team 
2021), the maps generated by MaxEnt were reclassified us- 
ing the 10-percentile training presence cloglog threshold, 
assuming an error rate of 10% in the training records, and 
rescaled between o and 1. Finally, all probability maps for 
the different habitat classes were combined by assigning 
each grid cell to the class with the highest probability. To 
remove some mottling for a smooth result, in R, the reso- 
lution of the map was resampled to a coarser resolution of 
50 X 50 cm per pixel using the nearest-neighbour approach 
with the function resample() from the raster package (Hy- 
MANS 2020). Based on the land cover mosaic map, distance 
maps quantifying the distance of every grid cell to all seven 
habitat categories were created (HIJMANS 2020). 


Principal Component Analysis 
and further statistical analysis 


In total, information on 32 habitat factors was collected 
for 947 encounters from 2018 to 2020. Dimensionality was 
reduced by performing a Principal Component Analy- 
sis (PCA). For this purpose, the data was sorted into four 
groups, including different habitat features comprised of 
functionally similar variables (weather, distance, microcli- 
mate, and microhabitat structure). Weather data is com- 
prised of the thirteen variables, current air temperature 
(in °C) 5 cm and 2 m above the ground, relative humid- 
ity (in %), amount (in mm) and duration of precipitation 
(in hours), solar irradiance (in joule), hours of sunshine, 
and maximum, minimum, and mean wind speed (in m/s) 
given for every ten minutes. Distances to habitat struc- 
tures are comprised of six variables, including shortest 
distances to sand, grass, heath, scrub, blackberry bushes, 
and trees. Microclimate is comprised of the five variables, 
cloud cover (in eighths), substrate temperature (in °C), air 
temperature (in °C), humidity (in %), and brightness of 
light at the encounter spot (sun = 1, half-shade = 2, shade = 
3). Microhabitat structure contains eight variables, sloping 
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(present = 1, not present = 0), open patches of sand, soil 
and/or gravel (present = 1, not present = o), grass (present 
= 1, not present = o), leaf litter (present =1, not present = 
o), low vegetation of grass and scrub (from 50 to 500 cm 
in height (present = 1, not present = o), trees taller than 
500 cm (present = 1, not present = o), substrate type (sand, 
gravel, soil and combinations), and the type of location 
(substrate = 1, grass = 2, leaf litter = 3, wood = 4, moss = 
5, heath = 6, blackberry = 7, or scrub = 8). The PCAs were 
conducted for all four groups separately and only Princi- 
pal Components (PCs) with eigenvalues > 1 were used for 
further analysis. The PCAs of the first two groups were 
performed with the function princomp() of the R package 
stats (R Core Team, 2021). The second two groups includ- 
ed also categorical variables and therefore, the function 
dudi.mix() of the package ade4 (DRAY & DUFOUR 2007) 
was used. 

Generalized Linear Models (GLMs), computed with the 
function glm() in R’s stats package (R Core Team 2021), 
were applied to test if the distance to habitat structures de- 
pended on time, weather, microhabitat structures, or mi- 
croclimatic conditions. The three different time variables, 
time in minutes since first encounter in 2018, Julian day, 
and calendar week were included to cover potential sea- 
sonal shifts. When a significant p-value was found (p < 
0.05), the test was repeated with only the significant terms 
and then the R? was calculated to examine how well the 
regression model fitted the observed data. The same ap- 
proach was applied to all microhabitat structure PCs, to 
investigate if microhabitat choice depended on sampling 
time or/and microclimatic conditions. 

To detect diel fluctuations in habitat use, density estima- 
tion plots using the function sm.density.compare() from 
the R package sm (BOWMAN & AZZALINI 2018) were con- 
ducted to show up the density of all distance PCs, micro- 
climate PCs, and microhabitat structure PCs for the three 
time-intervals, 8:00-11:00 h, 11:00-14:00 h, and 14:00- 
20:00 h. Additionally, pairwise permutation tests of equal- 
ity of the time intervals (morning/midday, morning/after- 
noon, midday/afternoon) were computed with 1,000 boot- 
straps to assess statistical significances (WILCOX 2010), and 
p-values were analyzed to assess if time intervals differed 
significantly from each other. 

Multidimensional hypervolumes (HVs) of the niche 
spaces were created using the function hypervolume_svm() 
of the R package hypervolume version 2.0.12 (BLONDER 
2019) to estimate niche differentiation between sexes (fe- 
male, male, subadult). This algorithm uses one-class sup- 
port vector machines (SVM) and is implemented with a 
radial basis function (RBF). The two parameters of the RBF 
were left at the default setting (y = 0.5 and v = 0.01), which 
is considered sufficient for such calculations (BLONDER 
2019). y determines an upper bound on the fraction of 
training errors and a lower bound of the fraction of sup- 
port vectors (lower values result in tighter wrapping of the 
shape to the data), and v defines the inverse radius of in- 
fluence of a single point (low values yield large influenc- 
es and smooth, less complex wraps around the data, and 
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high values yield small influences, make tighter but poten- 
tially noisier wraps around the data). The function hyper- 
volume_set() was used to visualize the niche of each sex 
class. Afterwards, the function hypervolume_overlap_sta- 
tistics() was used to compute pairwise overlap statistics for 
every comparison (females with males, females with sub- 
adults, and males with subadults). The function hypervol- 
ume_overlap_statistics() computes a set of four metrics: 
the Jaccard similarity that calculates the index of volume of 
intersection of 1 and 2 divided by the volume of union of 1 
and 2, the Sorensen similarity that calculates twice the vol- 
ume of intersection of 1 and 2 divided by volume of 1 plus 
the volume of 2, the unique fraction 1 that is the volume 
of unique component of 1 divided by volume of 1, and the 
unique fraction 2 that is the volume of unique component 
of 2 divided by volume of 2. To describe the influence of the 
PC axis on the HVs, the function hypervolume_variable_ 
importance() was used across 100 replicates for each com- 
parison, and the mean variable importance as well as the 
standard deviation were calculated. To investigate if sam- 
pling time (morning, midday, afternoon) influenced the 
encounter probability of a specific sex class, a Chi-squared 
test with the function chisq.test() available in the R pack- 
age stats (R Core Team 2021) was conducted. 

From 2018 to 2020, we obtained a total of 1095 Sand Liz- 
ard observations of which 947 had complete sets of meta- 
data and were used in the analysis. These were 429 adult 
females, 295 adult males, and 223 subadullts, i.e., immature 
juveniles and yearlings. We expected a low degree of pseu- 
doreplication as initial population size estimates based on 
photographic identification revealed very low re-encoun- 
ter rates (V. F. CLEMENT unpubl. data). 


Results 
Creating a categorical habitat map using MaxEnt 


Average AUC. „Values (Table 1) ranged from 0.655 for 
trees to 0.949 for heath. The AUC. of trees is compar- 
atively low, but still acceptable, as trees (predominantly 
Betulus pendula, Fagus sylvatica, Prunus serotina) cover a 
huge proportion of the study area and all other model per- 
formances exhibit good results, matching our experiences 
from the study site. 

The UTM tiles from 2019 turned out to be a very ac- 
curate representation of the overall habitat composition of 
the Dellbriicker Heide, compared to our on-site experience 
in 2020. The categorical map of the Dellbriicker Heide cov- 
ers an area of 36.8 ha, which represents almost the whole 
study area (40 ha) (Biindnis Heideterrasse e. V. 2021). Trees 
account for the most prominent vegetation class with over 
140,000 m? (Fig. 2). This is not surprising, because wood- 
lands frame all areas of lower vegetation, which mainly 
consist of grass and scrub while blackberry and heath are 
vegetation classes occurring at lower frequencies. Heath, 
the least frequent form of land cover, is only patchily dis- 
tributed, and blackberry occurs in a more scattered fash- 
ion across the entire area. Sandy patches are predominant 
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Table 1. Results of habitat classification calculated with MaxEnt. 


Class Training Averaged Threshold Band 1 Band 2 Band 3 
samples AUC,,,, + sd contribution contribution contribution 
Sand 1323 0.925 + 0.004 0.485 35.519 3.730 60.751 
Grass 4521 0.761 + 0.005 0.419 85.449 5.226 9.325 
Heath 214 0.949 + 0.010 0.375 53.174 44.620 2.206 
Scrub 1836 0.729 + 0.010 0.416 71.724 26.604 1.673 
Blackberry 157 0.810 + 0.026 0.451 64.356 24.133 11.511 
Trees 5944 0.655 + 0.006 0.490 71.336 12.174 16.490 
Water 8238 0.738 + 0.004 0.443 77.596 0.981 21.423 


in the Dellbriicker Heide where trails pass through and 
on the beach of the pond. Smaller artefacts are present in 
the transition area between sand and water. Shallow water, 
where the ground is visible, is wrongly classified as trees 
and grass. White areas in the map, which are mainly limit- 
ed to the canopy of trees, indicate that land cover classifica- 
tion was unsuccessful here. These areas were not included 
in the coverage calculation of the classes and therefore are 
neglected in the graphic of the mosaic map (Fig. 2). 


160,000 
140,000 
120,000 
100,000 
80,000 
60,000 
40,000 
20,000 


o 


Principal Component Analysis (PCA), 
Generalized Linear Models (GLMs) 
and density estimations 


The PCA based on weather conditions produced six PCs 
with eigenvalues > 1 (Table 2). The PCA based on distanc- 
es to habitat features produced five PCs with eigenvalues 
> ı (Table 3). The distances to water, as is represented by a 
large lake, was neglected because it had a very coarse gra- 


Figure 2. Land cover mosaic map depicting the seven habitat categories and the area covered by each in square metres (m°) with a 


resolution of 1 pixel = 50 x 50 cm. 
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Table 2. Summary of the Principal Component Analyses of weather conditions, showing factor loadings for each predictor, eigenvalues, 
explained variance, and the cumulative sums of the explained variance. The abbreviation a. g. stands for above ground. 


Variable PCI PC2 PC3 PC4 PC5 PC6 
Current temperature 200 cm a. g. (°C) 0.35 0.70 0.53 0.20 0.21 0.03 
Current temperature 5 cm a. g. (°C) 0.55 0.63 0.39 0.31 0.19 0.01 
Humidity 200 cm a. g. (%) -0.49 -0.70 0.51 0.03 0.03 0.00 
Maximum temperature (°C) 0.53 0.49 0.66 -0.21 -0.09 -0.04 
Minimum temperature (°C) 0.50 0.51 0.66 -0.20 -0.10 -0.04 
Minimum temperature 5 cm a. g. (°C) 0.73 0.36 0.54 -0.11 -0.12 0.00 
Solar irradiance (joule) 0.96 -0.27 -0.04 0.01 0.01 0.00 
Hours of sunshine 0.83 -0.08 -0.01 -0.15 -0.06 0.01 
Hours of precipitation -0.37 -0.05 -0.12 0.33 -0.03 -0.86 
Amount of precipitation (mm) -0.14 0.04 -0.08 0.10 0.04 -0.37 
Maximum wind speed (m/s) -0.08 -0.05 -0.12 -0.70 0.67 -0.08 
Minimum wind speed (m/s) -0.38 0.04 -0.12 -0.56 0.58 -0.06 
Mean wind speed (m/s) -0.21 0.01 -0.15 -0.67 0.65 -0.07 
Eigenvalues 294.69 100.75 62.10 9.51 6.16 2.69 
Explained variance 61.68 21.09 13.00 1.99 1.29 0.56 
Cumulative sum of exp. var. 61.68 82.77 95.76 97.75 99.04 99.61 


Table 3. Summary of the Principal Component Analyses of dis- 
tance, showing factor loadings for each predictor, eigenvalues, 
explained variance, and the cumulative sums of the explained 
variance. 


Table 5. Summary of the Principal Component Analysis of mi- 
crohabitat structure, showing factor loadings for each predictor, 
eigenvalues, explained variance, and the cumulative sums of the 
explained variance. 


Variable PC1 PC2 PC3 PC4 PC5 
Distance to sand 0.99 0.02 0.14 0.03 0.01 
Distance to grass 0.59 -0.09 -0.32 -0.63 -0.37 
Distance to scrub -0.13 -0.70 0.35 -0.44 0.37 
Distance to trees -0.32 -0.80 0.46 0.13 -0.16 
Distance to blackberry 0.24 -0.52 -0.06 -0.71 0.16 
Distance to heath 0.46 -0.63 -0.60 0.17 0.04 
Eigenvalues 38.16 11.86 6.94 3.49 1.24 
Explained variance 61.42 19.10 11.17 5.61 2.00 
Cumulative sum of 

exp. var. 61.42 80.51 91.69 97.30 99.30 


Table 4. Summary of the Principal Component Analysis of micro- 
climate, showing factor loadings for each predictor, eigenvalues, 
explained variance, and the cumulative sum of the explained vari- 
ance. “ifr” = measured with an infrared thermometer 


Variable PC1 PC2 
Clouds 0.48 0.58 
Substrate temperature (ifr) -0.79 0.00 
Air temperature -0.71 0.51 
Rel. humidity 0.76 -0.26 
Light condition at spot 0.37 0.77 
Eigenvalues 2.07 1:25 
Explained variance 41.49 25.08 
Cumulative sum of exp. var. 41.49 66.57 
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Variable PC1 PC2 PC3 
Slope 0.21 0.36 -0.08 
Open patches -0.09 0.56 0.59 
Grass 0.56 0.05 -0.51 
Leaf litter 0.06 -0.45 0.33 
Low vegetation -0.79 0.06 -0.25 
Trees -0.06 -0.55 0.45 
Substrate -0.16 0.49 0.25 
Location of observation -0.82 -0.07 -0.25 
Eigenvalues 1.70 1.21 1.12 
Explained variance 21.21 15.11 13.95 
Cumulative sum of exp. var. 21.21 36.32 50.27 


dient and was not directly present in our recording areas 
and therefore would have distorted the PCA results. Small, 
temporary water bodies were not mapped. The PCA based 
on microclimatic conditions produced two PCs with ei- 
genvalues > 1 (Table 4), and the PCA based on microhabi- 
tat structures produced three PCs with eigenvalues > 1 (Ta- 
ble 5). 

GLMs for some of the combinations produced signif- 
icant p-values (p < 0.05) (Supplementary Tables S3, 84), 
but when repeating the models with only the significant 
R° values were always lower than 0.07 (Supplementary Ta- 
bles S5, S6), suggesting very low performances of the re- 
gressions. 


Temporal habitat shifts in sand lizards? 


Density estimation plots for all distance-related PCs de- 
pending on the time intervals morning, midday, and af- 
ternoon display a slightly lower density in the afternoon 
for PC1 (Fig. 3). Distance PC1 is mainly composed of posi- 
tive correlations with the distances to sand and grass. The 


pairwise permutation tests of time groups also produced 
highly significant results for the comparisons of morning 
with afternoon (p = 0.004) and midday with afternoon 


(p < 0.001) (Fig. 3). 
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Figure 3. Density estimation plots for microclimate, distance, and microhabitat structure per time interval (morning/midday/after- 
noon), with p-values of pairwise permutation tests of equality. * p < 0.05, ** p < 0.01 and *** p < 0.005. 
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Table 6. Results of overlap statistics of the hypervolumes for pairwise comparisons of females, males and subadults. 


Group Pairwise Jaccard Sörensen Unique Unique 
comparison similarity similarity fraction 1 fraction 2 
Distance Female/male 0.538 0.700 0.318 0.281 
Female/ subadult 0.482 0.651 0.264 0.417 
Male/ subadult 0.467 0.636 0.258 0.443 
Weather Female/male 0.619 0.764 0.275 0.192 
Female/subadult 0.249 0.398 0.738 0.171 
Male/ subadult 0.245 0.394 0.734 0.243 
Microclimate Female/male 0.765 0.867 0.168 0.095 
Female/subadult 0.810 0.895 0.146 0.060 
Male/subadult 0.738 0.849 0.156 0.145 
Microhabitat Female/male 0.769 0.869 0.147 0.114 
structure Female/subadult 0.687 0.814 0.272 0.077 
Male/subadult 0.627 0.770 0.299 0.144 


As for microclimate, the density plot of PC1 shows slight 
differences between the time groups morning and midday 
and midday and afternoon, respectively (Fig. 3). Microcli- 
mate PC1 is mainly correlated with decreases of the sub- 
strate and air temperatures at the exact spot of an encoun- 
ter and an increase in humidity. The corresponding per- 
mutation tests with p-values of 0.005 and 0.003 also show 
highly significant results (Fig. 3), proving that densities at 
midday differ from those in the morning and in the after- 
noon. 

As for microhabitat structure, we only found differences 
depending on the time interval for PC2, showing a low- 
er density in the morning, whereas the density is slightly 
higher at midday and in the afternoon (Fig. 3). Microhab- 
itat structure PC2 contains mainly a positive correlation 
to open patches and a negative correlation to trees. Per- 
mutation tests reflect this (morning vs. midday, p = 0.024; 
morning vs. afternoon (p = 0.005; Fig. 3). Pairwise permu- 
tation tests of all other PCs were not significant. 


Hypervolumes and density estimation 


The overlap statistics of the hypervolumes (HVs) for pair- 
wise comparisons of females, males and subadults are pro- 
vided in Table 6. For distances, both similarity indices of all 
pairwise comparisons are highly similar. Additionally, both 
unique fractions of all pairwise comparisons of the HVs 
are low (< 32.0%), with only subadults showing unique 
fractions > 40.0% as compared to females and males. 

For weather, both similarity indices for pairwise com- 
parisons of HVs of females and males are very high, with a 
similarity of > 60.0%. The comparisons between the HVs 
of females or males and those of subadults reveal values 
< 40.0% for both similarity indices. The unique fractions 
of the HVs of females and males compared to those of sub- 
adults are also very high at values of 73.8% for females and 
73.4% for males. The unique fractions of HV of males com- 
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pared to those of females are quite low in contrast, reach- 
ing not more than 27.5% for females and 19.2% for males. 

For microclimate, both similarity indices for all pair- 
wise comparisons are very high, with a similarity of > 70%. 
Furthermore, unique fractions of the HVs of all sexes are 
very low, reaching not more than 16.8%. 

For microhabitat structures, both similarity indices for 
all pairwise comparisons are quite high with a similarity of 
> 68.5%. Furthermore, the niche of subadults does not dif- 
fer much from that of females (unique fraction of the HV 
of subadults only 7.0%). 

Density estimation plots for the sexes reveal an extreme- 
ly high density of data points within a very small range es- 
pecially for distance, which differ not differ between spe- 
cific sex classes (Fig. 4). Permutation tests reflect this by 
showing only two significant differences for distance PC2 
(p = 0.02) and PCs (p = 0.002) between males and females. 
Distance PC2 mainly contains negative correlations with 
distances to scrub, trees and heath. 

As for the factor weather, the density is quite broad. Per- 
mutation tests produced significant p-values for pairwise 
comparisons of males with subadults and females with 
subadults for weather PC1 and PC2. Weather PC1 con- 
tains mainly positive correlations with solar irradiance and 
hours of sunshine and minimum temperature 5 cm above 
the ground (Table 2). Weather PC2 includes a negative cor- 
relation with humidity and a positive correlation with tem- 
perature. Since both PCs share positive correlations with 
temperature, subadults appear to prefer different climatic 
conditions than do males and females. Apart from that, 
density estimation plots of weather PC1 suggest a bimod- 
al distribution of subadults. At weather PC4, which con- 
tains mainly minimum and maximum wind speeds (Ta- 
ble 2), the density plot shows small differences between 
sex classes (Fig. 4). Permutation tests are significant for the 
comparison of males with subadults (p = 0.01) and sug- 
gest that females share densities at the intersection of both 
other groups. Since weather PC6 only represents 0.5% of 
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Figure 4. Density estima- 
tion plots for weather, 
microclimate, distance, 
and microhabitat struc- 
ture per sex group with 
p-values of pairwise per- 
mutation tests of equality. 
* p < 0.05, ** p < 0.01 and 
*** p < 0.005. 
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the total variance of the weather variables the significant 
p-values were neglected. 

For microclimate PCi, mainly including the variable 
light conditions at spots, and microclimate PC2, mainly 
containing light conditions at the collection spots, cloud 
cover, and air temperatures, density estimation plots indi- 
cate that females have a lower density than males and sub- 
adults (Fig. 4). Permutation tests for microclimate do not 
result in any significant difference in densities between the 
sex groups at all. 

For microhabitat structure, density estimation plots 
show high similarities between all sexes for all PCs (Fig. 4). 
Permutation tests reveal only significant differences in mi- 
crostructure PC3 (p = 0.045) for males vs. females. Micro- 
structure PC3 mainly contains a positive correlation with 
open patches and trees and a negative correlation with 
grass coverage. 

Density estimation plots for the weather conditions on 
the sampling days (Fig. 5) illustrate that a broad range of 
different values are represented for every weather variable. 
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Permutation tests reflect this in the shape of significant dif- 
ferences for all weather variables (maximum temperature 
p = 0.000, relative humidity p = 0.000, mean wind speed 
= 0.03, hours of sunshine p = 0.000, hours of precipitation 
p = 0.000). 

The Chi-squared test for distribution, applied to detect 
if sampling time had influenced the encounter probability 
of a specific sex group, was not significant (p = 0.432, df = 
4, Chi-squared = 3.810). Therefore, sex groups and time in- 
tervals are statistically independent. 


Discussion 


Our results of potential seasonal habitat shifts of L. agilis in 
the Dellbriicker Heide, which is part of the distributional 
core area of the Sand Lizard, were compared with studies 
assessing populations from the species’ distributional range 
edges (JACKSON 1978, HOUSE & SPELLERBERG 1983, AMAT 
et al. 2003, NEMES et al. 2006, CEIRANS 2007, GROZDANOV 
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Figure 5. Density estimation plots for maximum temperature (°C), relative humidity (%), mean wind speed (m/s), hours of sunshine 
and precipitation (taken for every ten minutes), compared for all sampling days with p-values of permutation tests of equality. 


* p < 0.05 and ** p < 0.01. 
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et al. 2014, WOODFINE et al. 2017, PRIETO-RAMIREZ et al. 
2018). These studies suggest stenopotency for the assessed 
populations, likely depending on geographic scale and 
habitat type/availability. In contrast, we have proved here 
the eurypotency of L. agilis in a part of its distributional 
core by investigating habitat composition and factors with 
potential influence on habitat use. Despite our large sam- 
ple size and broad temporal coverage we found no support 
for the hypothesis that distances and microhabitat choice 
of Sand Lizards in the Dellbriicker Heide shift spatio-tem- 
porally according to daytime or season, which means that 
time had no influence on the movement of lizards in their 
habitat at all. Weather as well as microclimatic conditions 
apparently did not influence distances to habitat structures 
or the choice for structures in their microhabitat either. 
This indicates that the baseline conditions of the weather 
at the time of data acquisition (see Supplementary Fig. S1) 
in the Dellbriicker Heide did not induce lizards to engage 
in a specific behavioural pattern. We even found that the 
baseline conditions of the weather on sampling days were 
very representative and did not differ from random sam- 
pling. A high similarity in niche expression and niche sizes 
of lizards of all sex groups (male, female, subadult) exists 
for all categories (distance, weather, microclimate, micro- 
habitat structure), which suggests that lizards do not sig- 
nificantly differ as to their habitat use and we did not find 
significant differences in weather or microclimatic prefer- 
ences according to sex or size class, contrary to the results 
reported from distribution edges. Furthermore, sampling 
time did not affect our encountering of specific ontogenet- 
ic stages or sex groups. 


Influence of habitat factors on habitat use 


Our GLMs suggest that the Sand Lizards’ movement with- 
in their habitats is independent of the time, weather, struc- 
tures in their microhabitat, or microclimatic conditions. 
Density estimation plots and pairwise permutation tests 
showed that these lizards maintain a higher distance to 
sand and grass in the morning and midday, which is when 
they stay in the proximity of their hiding places in higher 
and denser vegetation. We also found slight differences in 
the microclimate when comparing midday with morning 
and afternoon. Temperatures in the morning are still low 
and lizards first must warm up in the morning and midday 
to stimulate their metabolism (Hous et al. 1980, AMAT et 
al. 2003, YANG et al. 2015). In the afternoon, when temper- 
atures decrease, these lizards are already warmed up and 
more mobile, for which reason they will then be more of- 
ten found in open areas composed of grass and sand. SAINT 
GrRrons (1976), who delved into activity periods of L. agilis 
in France, demonstrated that individuals left their hiding 
places early in the morning and became inactive already in 
the early afternoon, which matches our results. Although 
ecophysiological tests are required to prove this hypoth- 
esis, we expect that this might be caused by biotic factors. 
For instance, the lizards are more likely to have full stom- 


achs in the afternoon due to successful foraging earlier in 
the day and hence have fewer reasons to be active. What 
also needs to be considered is that humidity also will de- 
crease over time on a sunny day. Hence, habitats in the af- 
ternoon may warm up to similar temperatures but are dry- 
er. When we set into relation the time intervals with the 
number of individuals encountered we found that in total 
23.6% more individuals were spotted in the morning than 
in the afternoon. This is plausible when considering that 
lizards are still sluggish in the morning and are easier to 
spot while they bask at the bases of shrubs and in transi- 
tion areas to higher vegetation. This also is reflected by the 
results of the microhabitat structure data where the density 
for some coefficients is lower in the morning and at mid- 
day at PC2, which mainly is composed of a positive corre- 
lation to open patches in the habitat and negative to trees. 
This suggests that Sand Lizards more often sit near trees in 
the morning and at midday, but are more often surrounded 
by open patches in the afternoon. Nonetheless, these re- 
sults are only slightly obvious and not as pronounced as the 
density estimation plots suggest. 


Sex- and size-specific niche segregation 


Habitat niche utilization of the sex groups exhibit high 
similarities with only slight deviations in subadults. This 
might be due to the overall habitat composition of the Dell- 
briicker Heide. The whole area of 40 ha in total is quite 
small and contains several open areas that are framed by 
small pieces of woodland. Sand Lizards prefer to stay in the 
transition areas between lower vegetation and bare spots to 
higher vegetation (SCHIEMENZ & GUNTHER 1994, BLANKE 
2010, GROSSE & SEYRING 2015) and the Dellbriicker Hei- 
de has many of those transition areas, enabling these liz- 
ards to use them more extensively. Our results indicate that 
subadults use to a large extent the habitats occupied also 
by adult males and females, which is not surprising due to 
their being young males and females, but they also partly 
use other habitats to a smaller extent. This deviation could 
possibly be due to cannibalistic pressures on hatchlings 
when they still are very small. The latter is known to hap- 
pen in Sand Lizards, and yearlings may experience com- 
petitive pressure from adults. Juveniles and small subadults 
could be potential prey (BOHME 1984, CORBETT & TAMA- 
RIND 1979, PRIETO-RAMIREZ et al. 2018), and would be in- 
ferior opponents in interference competition (e.g., aggres- 
sion towards smaller individuals) as long as they are not 
large enough to assert themselves (DELANEY & WARNER 
2017, PRIETO-RAMIREZ et al. 2018). Another reason could 
be that juveniles, due to their smaller body sizes, are able 
to also use less expansive vegetation for hiding. Permuta- 
tion tests also show small differences in distances between 
males and females, which could be due to the breeding 
season, when males are known to move around more to 
find mating partners. This phase starts after hibernation in 
early spring and continues until July (GLANDT & BISCHOFF 
1988, BLANKE & FEARNLEY 2015). Female home ranges are 
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cited as often being much smaller than those of males, be- 
cause they might tend to stay nearer to possible nesting 
sites (OLSSON 1988, EDGAR & BIRD 2006). This is also rep- 
resented in permutation tests for microhabitat structure 
that reveal a broader variation in many habitat structures 
in males than in females. Such a pattern can be expected, as 
mate guarding and home range defence is well known from 
the males of this species (OLSON 1993, OLSON et al. 1996). 

According to the weather HVs, males and females are 
highly similar in niches with small unique fractions, show- 
ing that weather-dependent niches of males and females do 
not differ whereas subadults differ slightly in their weather 
preferences from those of males and females. This is also 
reflected in the density estimation plot of weather relative 
to sex group, according to which subadults do not only de- 
viate partly from adults but also show a bimodal distribu- 
tion for PC1 and PC2 containing solar irradiance, hours 
of sunshine, temperature, and humidity. These parame- 
ters are very important for thermoregulation in Sand Liz- 
ards (JACKSON 1978, DENT & SPELLERBERG 1987, AMAT et 
al. 2003). Younger individuals have smaller volumes and 
therefore heat up faster, which means that they face a high- 
er risk of overheating, but cannot on the other hand keep 
body temperatures constant for longer periods, forcing 
them to alternate between basking and shaded spots more 
frequently than adults, which heat up more slowly, but due 
to their larger volumes can keep temperatures within the 
preferred range for longer periods, therefore have to bask 
less often, and can exploit colder areas for longer periods 
of time (Saconas et al. 2013). BÖHME (1984) stated that the 
preferred body temperature of subadults was slightly lower 
than that of adults. The bimodal distribution in subadults 
might be an artefact due to the category subadults not only 
comprising small juveniles but also yearlings, which are al- 
ready quite large and share more similarities with the adult 
males and females. The deviation in subadults in density 
estimation plots is consistent with the HVs of weather and 
distance. Overall, Sand Lizards appear to occupy large cli- 
matic niches for specific habitat factors independent of 
their sex. This is evident from the broader spread of the 
densities for weather in PC1 in the respective estimation 
plots. The slight differences in size and sex in PCi in the 
weather density plot can be explained by the fact that ju- 
veniles due to their smaller sizes thermoregulate more fre- 
quently. Sand Lizards in the Dellbriicker Heide are active 
not only in similar weather conditions, but also do not 
show any specializations for specific microclimatic condi- 
tions, as is suggested not only by the respective HVs but 
also permutation tests for equality. 

Furthermore, when considering the choice of micro- 
habitat structures by these lizards, HVs of sexes again come 
up with high similarities in all groups, and the fact that 
subadults have nearly no unique fractions compared to fe- 
males shows that subadults use the same but fewer micro- 
habitat structures and thus have a smaller niche size com- 
pared with both males or females, proving that they avoid 
certain areas in which they are out-competed by adults or 
are at risk of filial cannibalism (CORBETT & TAMARIND 
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1979, BOHME 1984, PRIETO-RAMIREZ et al. 2018). Density 
estimation plots of microhabitat structure are consistent 
with the results of the HVs in that they indicate no sig- 
nificant differences in densities for subadults. Density plots 
overall illustrate a wide range in densities, evidencing that 
these lizards make use of a high diversity of microhabitat 
structures. Our study results also prove that the lack of dif- 
ferences in habitat use in the Dellbriicker Heide could not 
be caused by a differentiation between sex groups along 
diel timelines by showing that encounter times had no in- 
fluence on the spotting of a specific sex group. 


Core versus periphery of distribution 


All statistic evaluations from the Dellbriicker Heide pop- 
ulation of the Sand Lizard, which is part of the distribu- 
tional core area of this species, provide completely differ- 
ent results to those reported from the edges/periphery of 
its distributional range, from which differences were re- 
ported to be very pronounced and influenced by various 
habitat characteristics (JACKSON 1978, HOUSE & SPELLER- 
BERG 1983, AMAT et al. 2003, NEMES et al. 2006, CEIRANS 
2007, GROZDANOV et al. 2014, WOODFINE et al. 2017, PRIE- 
TO-RAMIREZ et al. 2018). Our results show only slight dif- 
ferences in habitat use by the Sand Lizards, which are at- 
tributed only to sexes and ontogenetic stage, meaning that 
these differences are shaped by characteristics of the liz- 
ards themselves but not by external circumstances. HENLE 
et al. (2017) suggested that populations at the range edges 
of the Sand Lizard expressed a higher degree of speciali- 
zation, had a lower diversity, and were more sensitive to 
habitat fragmentation compared to those at the core. This 
is reasonable when assuming that in the core region, living 
conditions for these lizards are optimal due to the com- 
bination of a wider range of suitable microhabitats (HEN- 
LE et al. 2017, PRIETO-RAMIREZ et al. 2020). This in turn 
makes the lizards become habitat generalists in their distri- 
butional core area, compared to the periphery where they 
have to adopt a higher degree of specialization to compen- 
sate for altogether poorer habitat conditions (BOHME 1978, 
SCHAEFER 2012, BOHME & RODDER 2014, PRIETO-RAMI- 
REZ et al. 2018). Our results prove a high range of tolerance 
towards many of the evaluated habitat factors of L. agilis 
in this population. Compared to populations studied in 
peripheral distribution ranges where lizards tended to be 
bound to specific habitats to an extreme extent, this con- 
stitutes an enormous contrast and shows this species’ great 
capacity of behavioral adjustment to environmental con- 
ditions. In the absence of evidence for local adaptation or 
pronounced acclimatization the more parsimonious hy- 
pothesis of more limiting factors at the range edges com- 
pared to the core should be preferred. The latter is known 
as the KUHNELT (1943) principle, which operates at popu- 
lation level and is well applicable to L. agilis. The habitat 
mosaic of the Dellbriicker Heide is a perfect representa- 
tion of the distributional core area of the Sand Lizard, pro- 
viding optimal living conditions and a very pronounced 
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landscape heterogeneity, which is proven by none of the 
evaluated habitat factor having an influence on habitat use. 
Therefore, proper management of those core areas with the 
aim of preserving habitat variability, as well as re-establish- 
ing connectivity of fragmented habitats especially in the 
distribution edges, are very important aspects for the con- 
servation of the entire Sand Lizard population. 
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